from django.http import HttpResponse, Http404, HttpResponseRedirect
from django.template import TemplateDoesNotExist, RequestContext
from django.shortcuts import render_to_response
from django import forms
from settings import MEDIA_ROOT
from numpy import asarray, zeros,matrix
import networkx as nx
from hcluster import *
import pickle
import time

class netdecform(forms.Form):
    maxinisize=forms.CharField(label="Max size of initial modules", required=False)
    minsize=forms.CharField(label="Minimal size of modules", required=False)
    nm=forms.CharField(label="final number of modules", required=False)
    file  = forms.FileField()
#    pa=forms.ChoiceField(label="select pathways", widget=forms.SelectMultiple, choices=p_options)

def netdec(request):
    if request.method == 'POST':
        f = netdecform(request.POST,request.FILES)
        if not f.is_valid():
           err="sth wrong"
           return render_to_response('netdec.html', locals(),context_instance=RequestContext(request)) #do sth else 
        else:
           if f.cleaned_data["maxinisize"]:
             threshold=int(f.cleaned_data["maxinisize"])
           else:
             threshold=30
           if f.cleaned_data["minsize"]:
             s=int(f.cleaned_data["minsize"])
           else:
             s=0
           if f.cleaned_data["nm"]:
             nmod=int(f.cleaned_data["nm"])
           else:
             nmod=1
           if f.cleaned_data["file"]:
             net=request.FILES['file']
             line=net.readline()
             if '\t' in line:
                G=nx.read_edgelist(net,delimiter='\t')
                a=len(G)
             elif ',' in line:
                G=nx.read_edgelist(net,delimiter=',')
                a=len(G)
             else:   
                G=nx.read_edgelist(net) #by default delimited by space
                b=len(G)
#             G=nx.read_edgelist(net,delimiter='\t',data=(('lable',str),)) #tab delimited and with other information
             UG=G.to_undirected()
             compsub=nx.connected_component_subgraphs (UG)
             size=[len(wc) for wc in compsub]
             m=max(size)
             igwc=size.index(m)
             gwc=compsub[igwc]
             mnodes=treemod(gwc,threshold)
             if s!=0: #use first optimization to merge small modules
               mnodes=firstoptmod(gwc,mnodes,s)
             mdict,Qlist=secondopt1re(gwc,mnodes,nmod)
             f3=open(MEDIA_ROOT+'/decres.txt','w')
             sn='number'
             sl='min size'
             sq='modularity'
             for items in Qlist:
               sn=sn+'\t'+str(items[0])
               sl=sl+'\t'+str(items[1])
               sq=sq+'\t'+str(items[2])[:7]
             f3.write(sn+'\n')
             f3.write(sl+'\n')
             f3.write(sq+'\n')
             for k,v in mdict.iteritems():
               f3.write('%s\t%s\n' % (k, v))
             f3.close()  
             res=open(MEDIA_ROOT+'/decres.txt',"rb").read()
             response=HttpResponse(res, mimetype="text/plain")
             response['Content-Disposition']='attachment; filename=realist.txt'
             return response
#             return render_to_response('krea.html', {'reaction_list': objs})
           else:
               err="please choose the network file!"
               return render_to_response('netdec.html', locals(),context_instance=RequestContext(request))
    else:
        f = netdecform(request.POST)
        err="please choose the network file!"
        return render_to_response('netdec.html', locals(),context_instance=RequestContext(request))

def treemod(G,threshold):
  gwcnodes=G.nodes()
  gwcnodes.sort()
  ngwc=len(gwcnodes)
  length=nx.all_pairs_shortest_path_length(G)
  aa=[[l[k] for k in sorted(l.keys())] for l in [length[k] for k in sorted(length.keys())]]
  x=asarray(aa)
  m=matrix(x)   # the distance matrix
  plength=squareform(m) # change the distance matrix into the condensed format
  Z=linkage(plength,method='average')
#  Z=linkage(plength,method='complete') #not affect the result
  tree=totree(Z)
# clustering of the tree with threshold=50 and put it into file
  mycurNode = [None] * (2 * ngwc)
  lvisit = zeros((2 * ngwc,), dtype=bool)
  rvisit = zeros((2 * ngwc,), dtype=bool)
  mycurNode[0] = tree
  myk = 0
  mnodes=[]
  while myk >= 0:
        mynd = mycurNode[myk]
        myndid = mynd.id
        if not mynd.isLeaf():
            if not lvisit[myndid]:
                lvisit[myndid] = True
                lc=mynd.left.count #number of leaf nodes below this node
                if lc<threshold:
                    k = myk+1 #level from the current node, the following part is modified from original preOrder function
#                    mid=mid+1
                    mycurNode[k]=mynd.left
                    m=[]
                    while k > myk:
                        nd = mycurNode[k]
                        ndid = nd.id
                        if nd.isLeaf():
                            m.append(gwcnodes[ndid])
#                            f.write(str(gwcnodes[ndid])+'\t'+str(igwc)+'\t'+str(mid)+'\n')
                            k = k - 1
                        else:
                            if not lvisit[ndid]:
                                mycurNode[k + 1] = nd.left
                                lvisit[ndid] = True
                                k = k + 1
                            elif not rvisit[ndid]:
                                mycurNode[k + 1] = nd.right
                                rvisit[ndid] = True
                                k = k + 1
                            # If we've visited the left and right of this non-leaf
                            # node already, go up in the tree.
                            else:
                                k = k - 1
                    mnodes.append(m)
                    myk=myk-1
                else:
                    mycurNode[myk + 1] = mynd.left
                myk = myk + 1
            elif not rvisit[myndid]:
                rvisit[myndid] = True
                rc=mynd.right.count
                if rc<threshold:
                    k = myk+1 #level from the current node, the following part is modified from original preOrder function
#                    mid=mid+1
                    mycurNode[k]=mynd.right
                    m=[]
#                    curNode[0]=mynd.right
                    while k > myk:
                        nd = mycurNode[k]
                        ndid = nd.id
                        if nd.isLeaf():
                            m.append(gwcnodes[ndid])
#                            f.write(str(gwcnodes[ndid])+'\t'+str(igwc)+'\t'+str(mid)+'\n')
                            k = k - 1
                        else:
                            if not lvisit[ndid]:
                                mycurNode[k + 1] = nd.left
                                lvisit[ndid] = True
                                k = k + 1
                            elif not rvisit[ndid]:
                                mycurNode[k + 1] = nd.right
                                rvisit[ndid] = True
                                k = k + 1
                            # If we've visited the left and right of this non-leaf
                            # node already, go up in the tree.
                            else:
                                k = k - 1
                    mnodes.append(m)
                    myk=myk-1
                else:
                    mycurNode[myk + 1] = mynd.right
                myk = myk + 1
                # If we've visited the left and right of this non-leaf
                # node already, go up in the tree.
            else:
                myk = myk - 1
  return mnodes

def modularity(G,nodelist): #better
    L=float(G.number_of_edges())
    Q=0.0
    n=len(nodelist)
    ls=[len(nx.edge_boundary(G,nodelist[i],nodelist[i])) for i in range(n)]
    ds2=[(sum(G.degree(nodelist[i]).values()))**2 for i in range(n)]
    Q=0.5*sum(ls)/L-sum(ds2)/4.0/L/L
    return Q

def firstoptmod(G,nodelist,threshold): #better result than just considering boundary
#combine the small cluster with less than 10 nodes into bigger clusters, based on modularity
  nodelist.sort(cmp=bylength) 
  L=float(G.number_of_edges())
  while len(nodelist[0])<threshold:
        mc=-0.1 #the max inter-module links
        # check all the clusters which connect with cluster i
        for j in xrange(1,len(nodelist)):  
            k=len(nx.edge_boundary(G,nodelist[0],nodelist[j]))
            if k>0:
               Mchange=k/L-sum(G.degree(nodelist[0]).values())*sum(G.degree(nodelist[j]).values())/2.0/L/L
               if Mchange>mc:
                    comb=j
                    mc=Mchange
        nodelist[comb]=nodelist[0]+nodelist[comb]
        nodelist.remove(nodelist[0])
        nodelist.sort(cmp=bylength)
  return nodelist
  
def bylength(a,b):
 	return len(a)-len(b)

def secondopt1re(G,tpndlist,num): #the best
  #using a reduced network graph to reduce the time for the loops
  tpndlist.sort(cmp=bylength)  # order the clusters by size, will make the program run fast
  Qlist=[]
  Q=modularity(G,tpndlist)
  Qlist.append((len(tpndlist),len(tpndlist[0]),Q)) #saving the number of modules and the modularity
  L=float(G.number_of_edges())
  Gre=nx.blockmodel(G,tpndlist)
  '''    
  a=[]
  for l in tpndlist:
    a=a+[(x,tpndlist.index(l)) for x in l]
  mdict=dict(a)
  '''
  mdict={}  #saving a series of partitions as dictionary
  for mi in range(len(tpndlist)):
    for node in tpndlist[mi]:
      mdict[node]=str(mi+1)
  while len(tpndlist)>num:
    mc=-0.1
    for e in Gre.edges():
               k=len(nx.edge_boundary(G,tpndlist[e[0]],tpndlist[e[1]]))
               Mchange=k/L-sum(G.degree(tpndlist[e[0]]).values())*sum(G.degree(tpndlist[e[1]]).values())/2.0/L/L
#                    Mchange=modularity_change(G,tpndlist,i,j)
                # find the proper cluster with largest modularity change and smallest size
               if Mchange>mc:
                        comb=e
                        mc=Mchange
    tpndlist[comb[1]]=tpndlist[comb[0]]+tpndlist[comb[1]]
    tpndlist.remove(tpndlist[comb[0]])
    tpndlist.sort(cmp=bylength)  # order the clusters by size
    Gre=nx.blockmodel(G,tpndlist)
    '''    
    a=[]
    for l in tpndlist:
      a=a+[(x,tpndlist.index(l)) for x in l]
    mdict=dict(a)
    '''
    for mi in xrange(len(tpndlist)):
      for node in tpndlist[mi]:
        mdict[node]='%s\t%s' % (mdict[node],mi+1)
    Q=Q+mc    
    Qlist.append((len(tpndlist),len(tpndlist[0]),Q))
  return mdict,Qlist

